To quantitatively examine the efficacy of vegetation restoration in drylands globally.
Study-level viz to document patterns in exclusions primarily and the relatie frequenices, at the study level, of major categories of evidence.
#study data####
library(tidyverse)
studies <- read_csv("data/studies.csv")
studies
## # A tibble: 278 x 18
## ID title technique data region exclude rationale observations
## <dbl> <chr> <chr> <chr> <chr> <chr> <chr> <chr>
## 1 152 Shor… seeding,… expe… Africa no <NA> <NA>
## 2 180 Rest… chemical… App.… Africa no <NA> <NA>
## 3 229 Infl… soil see… fiel… Africa no <NA> <NA>
## 4 230 Acti… planting fiel… Africa no <NA> <NA>
## 5 255 The … grazing … fiel… Africa no <NA> <NA>
## 6 262 Reve… seeding,… eper… Africa no <NA> <NA>
## 7 263 The … phytogen… fiel… Africa no <NA> <NA>
## 8 264 Eval… seeding,… fiel… Africa no <NA> <NA>
## 9 271 Patc… natural … fiel… Africa no <NA> <NA>
## 10 4 Fact… natural … App.… Africa no <NA> <NA>
## # ... with 268 more rows, and 10 more variables: disturbance <chr>,
## # system <chr>, goal <chr>, intervention <chr>, paradigm <chr>,
## # grazing <chr>, hypothesis <chr>, soil <chr>, benchmark <chr>,
## # notes <chr>
#quick look at rationale needed
exclusions <- studies %>%
filter(exclude == "yes")
#quick look at studies with paradigms
evidence <- studies %>%
filter(exclude == "no")
#library(skimr)
#skim(evidence)
#study-level viz#####
#exclusions
ggplot(exclusions, aes(rationale, fill = region)) +
geom_bar() +
coord_flip() +
labs(x = "rational for exclusion", y = "frequency") +
scale_fill_brewer(palette = "Paired")
ggplot(evidence, aes(disturbance, fill = paradigm)) +
geom_bar(na.rm = TRUE) +
coord_flip() +
scale_fill_brewer(palette = "Paired") +
labs(y = "frequency")
ggplot(evidence, aes(region, fill = paradigm)) +
geom_bar(na.rm = TRUE) +
coord_flip() +
scale_fill_brewer(palette = "Paired") +
labs(y = "frequency")
ggplot(evidence, aes(data, fill = paradigm)) +
geom_bar(na.rm = TRUE) +
coord_flip() +
scale_fill_brewer(palette = "Paired") +
labs(y = "frequency")
ggplot(evidence, aes(system, fill = paradigm)) +
geom_bar(na.rm = TRUE) +
coord_flip() +
scale_fill_brewer(palette = "Paired") +
labs(y = "frequency")
ggplot(evidence, aes(goal, fill = paradigm)) +
geom_bar(na.rm = TRUE) +
coord_flip() +
scale_fill_brewer(palette = "Paired") +
labs(x = "outcome", y = "frequency")
#step 1 models####
#paradigm
derived.evidence <- evidence %>%
group_by(technique, data, region, disturbance, goal, paradigm) %>% summarise(n = n())
#active-passive split
m <- glm(n~paradigm, family = poisson, derived.evidence)
anova(m, test="Chisq")
## Analysis of Deviance Table
##
## Model: poisson, link: log
##
## Response: n
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 167 9.9147
## paradigm 1 0.045115 166 9.8696 0.8318
#region
m1 <- glm(n~paradigm*region, family = poisson, derived.evidence)
#m1
#summary(m1)
anova(m1, test="Chisq")
## Analysis of Deviance Table
##
## Model: poisson, link: log
##
## Response: n
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 167 9.9147
## paradigm 1 0.045115 166 9.8696 0.8318
## region 6 0.301367 160 9.5682 0.9995
## paradigm:region 6 0.213627 154 9.3546 0.9998
#outcome
m2 <- glm(n~paradigm*goal, family = poisson, derived.evidence)
#m1
#summary(m1)
anova(m2, test="Chisq")
## Analysis of Deviance Table
##
## Model: poisson, link: log
##
## Response: n
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 167 9.9147
## paradigm 1 0.045115 166 9.8696 0.8318
## goal 6 0.240941 160 9.6287 0.9997
## paradigm:goal 4 0.301480 156 9.3272 0.9897
#even split between active and passive evidence by all key categories
A summary of sort process using PRISMA.
library(PRISMAstatement)
prisma(found = 1504,
found_other = 5,
no_dupes = 1039,
screened = 1039,
screen_exclusions = 861,
full_text = 178,
full_text_exclusions = 100,
qualitative = 100,
quantitative = 78,
width = 800, height = 800)
Check data and calculate necessary measures.
data <- read_csv("data/data.csv")
data <- data %>%
mutate(lrr = log(mean.t/mean.c), rii = ((mean.t-mean.c)/(mean.t + mean.c)), var.es = ((sd.t^2/n.t*mean.t^2) + (sd.c^2/n.c*mean.c^2)))
data
## # A tibble: 3,291 x 48
## study.ID ID publication.year data.year exp.year `exp.length (mo…
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 1.1 2018 NA NA NA
## 2 1 1.1 2018 NA NA NA
## 3 1 1.1 2018 NA NA NA
## 4 1 1.2 2018 NA NA NA
## 5 1 1.2 2018 NA NA NA
## 6 1 1.2 2018 NA NA NA
## 7 1 1.3 2018 NA NA NA
## 8 1 1.3 2018 NA NA NA
## 9 1 1.3 2018 NA NA NA
## 10 1 1.4 2018 NA NA NA
## # ... with 3,281 more rows, and 42 more variables: disturbance <chr>,
## # focus <chr>, technique <chr>, paradigm <chr>, hypothesis <chr>,
## # pathway <chr>, plant.species <chr>, target.plant <chr>,
## # measure.success <chr>, measured.factor <chr>, factor.levels <chr>,
## # treatment <chr>, control <chr>, unit <chr>, Nsites <dbl>, n.t <dbl>,
## # n.c <dbl>, ntotalsamples <dbl>, mean.t <dbl>, mean.c <dbl>,
## # sd.t <dbl>, sd.c <dbl>, se.t <dbl>, se.c <dbl>, `p-value` <dbl>,
## # df <dbl>, measure.dispersion <chr>, lat <dbl>, long <dbl>,
## # continent <chr>, country <chr>, ecosystem <chr>, `elevation
## # (m)` <dbl>, `MAP (mm)` <dbl>, aridity.index <dbl>,
## # `potential.evaporation (mm)` <dbl>, grazing <chr>, soil <chr>,
## # notes <chr>, lrr <dbl>, rii <dbl>, var.es <dbl>
#consider adding some other effect size measures and/or study-level data too
Explore summary level data of all data. Explore aggregation levels that support the most reasonable data structure and minimize non-independence issues.
#evidence map####
require(maps)
world<-map_data("world")
map<-ggplot() + geom_polygon(data=world, fill="gray50", aes(x=long, y=lat, group=group))
map + geom_point(data=data, aes(x=long, y=lat)) #render a literal map, i.e. evidence map, of where we study the niche in deserts globally
#aggregation####
data.simple <- data %>%
group_by(study.ID, paradigm, technique, measure.success) %>%
summarise(n = n(), mean.lrr = mean(lrr), mean.rii = mean(rii), mean.var = mean(var.es))
se <- function(x){
sd(x)/sqrt(length(x))
}
simple.data <- data %>% group_by(study.ID, paradigm) %>% summarise(mean.rii = mean(rii), error = se(rii))
#viz for aggregation####
ggplot(na.omit(data.simple), aes(technique, n, fill = paradigm)) +
geom_bar(stat = "identity") +
coord_flip() +
scale_fill_brewer(palette = "Paired")
ggplot(na.omit(data.simple), aes(measure.success, n, fill = paradigm)) +
geom_bar(stat = "identity") +
coord_flip() +
scale_fill_brewer(palette = "Paired")
Meta and conventional statistical models to explore relative efficacy.
library(plotrix) #for quick s.e. calculations sometimes needed for data tidy step
library(meta) #nice package for most meta-statistics
#assign model (typically a nice meta. function from one of several packages such as meta, metafor, or netmeta)
m <- metagen(mean.rii, error, studlab = study.ID, byvar = paradigm,data = simple.data) #fit generic meta-analysis to an object
m
## 95%-CI %W(fixed) %W(random) paradigm
## 1 0.1004 [-0.0416; 0.2424] 0.5 3.9 passive
## 2 NA 0.0 0.0 passive
## 3 -0.1868 [-0.2919; -0.0817] 1.0 4.2 active
## 4 NA 0.0 0.0 passive
## 5 NA 0.0 0.0 active
## 6 0.1973 [-0.1756; 0.5702] 0.1 2.0 active
## 9 0.1901 [ 0.1523; 0.2280] 7.5 4.6 passive
## 11 NA 0.0 0.0 active
## 12 NA 0.0 0.0 active
## 13 -0.0083 [-0.0913; 0.0747] 1.6 4.3 active
## 15 NA 0.0 0.0 passive
## 16 NA 0.0 0.0 active
## 19 NA 0.0 0.0 active
## 20 -0.2808 [-0.3514; -0.2102] 2.2 4.4 active
## 21 NA 0.0 0.0 active
## 23 NA 0.0 0.0 passive
## 28 NA 0.0 0.0 passive
## 29 -0.2971 [-0.3868; -0.2073] 1.3 4.3 active
## 30 -0.0326 [-0.1910; 0.1257] 0.4 3.7 active
## 31 NA 0.0 0.0 passive
## 32 NA 0.0 0.0 active
## 34 NA 0.0 0.0 passive
## 35 NA 0.0 0.0 passive
## 37 NA 0.0 0.0 passive
## 39 NA 0.0 0.0 passive
## 40 NA 0.0 0.0 active
## 41 NA 0.0 0.0 passive
## 42 NA 0.0 0.0 active
## 43 NA 0.0 0.0 passive
## 44 NA 0.0 0.0 active
## 45 NA 0.0 0.0 active
## 46 NA 0.0 0.0 passive
## 48 NA 0.0 0.0 active
## 50 NA 0.0 0.0 active
## 51 0.1280 [-0.0646; 0.3205] 0.3 3.4 active
## 52 NA 0.0 0.0 active
## 53 NA 0.0 0.0 active
## 54 NA 0.0 0.0 active
## 55 NA 0.0 0.0 active
## 58 NA 0.0 0.0 passive
## 59 NA 0.0 0.0 active
## 61 NA 0.0 0.0 active
## 63 NA 0.0 0.0 active
## 64 NA 0.0 0.0 active
## 65 NA 0.0 0.0 active
## 66 0.0762 [ 0.0222; 0.1302] 3.7 4.5 passive
## 68 -0.1772 [-0.2845; -0.0699] 0.9 4.2 active
## 69 0.5020 [ 0.3552; 0.6489] 0.5 3.8 passive
## 71 NA 0.0 0.0 passive
## 73 NA 0.0 0.0 passive
## 74 NA 0.0 0.0 active
## 76 NA 0.0 0.0 active
## 77 0.0613 [-0.8629; 0.9854] 0.0 0.5 active
## 78 NA 0.0 0.0 passive
## 80 NA 0.0 0.0 active
## 82 NA 0.0 0.0 passive
## 84 NA 0.0 0.0 active
## 85 NA 0.0 0.0 active
## 86 NA 0.0 0.0 active
## 87 -0.0138 [-0.0406; 0.0131] 14.8 4.6 active
## 88 -0.1297 [-0.1597; -0.0997] 11.9 4.6 active
## 89 NA 0.0 0.0 active
## 91 NA 0.0 0.0 passive
## 93 NA 0.0 0.0 passive
## 95 -0.2280 [-0.3094; -0.1466] 1.6 4.3 active
## 96 NA 0.0 0.0 passive
## 98 NA 0.0 0.0 active
## 100 NA 0.0 0.0 active
## 101 NA 0.0 0.0 passive
## 104 0.3069 [ 0.2619; 0.3519] 5.3 4.5 active
## 105 NA 0.0 0.0 passive
## 106 NA 0.0 0.0 active
## 107 NA 0.0 0.0 active
## 108 NA 0.0 0.0 passive
## 109 0.4675 [ 0.2613; 0.6737] 0.3 3.3 active
## 110 NA 0.0 0.0 passive
## 111 0.0748 [ 0.0406; 0.1091] 9.2 4.6 active
## 112 NA 0.0 0.0 active
## 114 NA 0.0 0.0 active
## 115 NA 0.0 0.0 passive
## 118 NA 0.0 0.0 active
## 119 NA 0.0 0.0 active
## 120 NA 0.0 0.0 passive
## 121 -0.0229 [-0.0432; -0.0026] 26.0 4.6 active
## 126 NA 0.0 0.0 passive
## 127 NA 0.0 0.0 active
## 128 NA 0.0 0.0 passive
## 129 NA 0.0 0.0 active
## 130 NA 0.0 0.0 passive
## 132 NA 0.0 0.0 active
## 133 NA 0.0 0.0 passive
## 134 NA 0.0 0.0 active
## 135 0.0859 [ 0.0317; 0.1401] 3.7 4.5 passive
## 139 NA 0.0 0.0 active
## 141 NA 0.0 0.0 passive
## 142 NA 0.0 0.0 active
## 143 NA 0.0 0.0 active
## 145 NA 0.0 0.0 active
## 146 NA 0.0 0.0 active
## 147 -0.2642 [-0.4402; -0.0882] 0.3 3.6 active
## 148 NA 0.0 0.0 active
## 149 NA 0.0 0.0 active
## 151 NA 0.0 0.0 active
## 152 NA 0.0 0.0 active
## 153 NA 0.0 0.0 active
## 154 NA 0.0 0.0 active
## 155 NA 0.0 0.0 active
## 157 NA 0.0 0.0 active
## 160 NA 0.0 0.0 active
## 161 NA 0.0 0.0 active
## 162 NA 0.0 0.0 passive
## 164 NA 0.0 0.0 active
## 166 NA 0.0 0.0 active
## 167 NA 0.0 0.0 active
## 168 NA 0.0 0.0 active
## 169 NA 0.0 0.0 active
## 170 NA 0.0 0.0 passive
## 171 NA 0.0 0.0 active
## 179 NA 0.0 0.0 active
## 180 NA 0.0 0.0 active
## 182 NA 0.0 0.0 active
## 188 NA 0.0 0.0 passive
## 193 NA 0.0 0.0 active
## 194 NA 0.0 0.0 active
## 195 NA 0.0 0.0 active
## 197 NA 0.0 0.0 active
## 202 NA 0.0 0.0 active
## 204 NA 0.0 0.0 active
## 206 NA 0.0 0.0 active
## 207 NA 0.0 0.0 passive
## 209 NA 0.0 0.0 passive
## 210 -0.3286 [-0.4648; -0.1924] 0.6 3.9 active
## 211 NA 0.0 0.0 active
## 212 NA 0.0 0.0 active
## 213 NA 0.0 0.0 active
## 214 NA 0.0 0.0 active
## 218 NA 0.0 0.0 active
## 219 NA 0.0 0.0 passive
## 220 NA 0.0 0.0 passive
## 221 NA 0.0 0.0 passive
## 223 NA 0.0 0.0 passive
## 229 NA 0.0 0.0 passive
## 230 NA 0.0 0.0 active
## 231 NA 0.0 0.0 passive
## 232 NA 0.0 0.0 active
## 233 NA 0.0 0.0 active
## 235 NA 0.0 0.0 active
## 238 NA 0.0 0.0 active
## 239 0.2968 [-0.1971; 0.7907] 0.0 1.4 passive
## 240 NA 0.0 0.0 active
## 241 NA 0.0 0.0 active
## 242 NA 0.0 0.0 active
## 243 NA 0.0 0.0 active
## 244 NA 0.0 0.0 active
## 245 NA 0.0 0.0 active
## 246 NA 0.0 0.0 active
## 247 -0.3379 [-0.3809; -0.2950] 5.8 4.5 passive
## 248 NA 0.0 0.0 active
## 250 NA 0.0 0.0 active
## 251 NA 0.0 0.0 active
## 252 NA 0.0 0.0 active
## 253 NA 0.0 0.0 active
## 254 NA 0.0 0.0 active
## 255 NA 0.0 0.0 passive
## 256 0.0430 [-0.0941; 0.1801] 0.6 3.9 passive
## 259 NA 0.0 0.0 passive
## 262 NA 0.0 0.0 active
## 263 NA 0.0 0.0 passive
## 264 NA 0.0 0.0 active
## 267 NA 0.0 0.0 passive
## 268 NA 0.0 0.0 active
## 269 NA 0.0 0.0 active
## 271 NA 0.0 0.0 passive
## 274 NA 0.0 0.0 active
## 275 NA 0.0 0.0 active
## 276 NA 0.0 0.0 active
## 277 NA 0.0 0.0 active
## 278 NA 0.0 0.0 active
##
## Number of studies combined: k = 26
##
## 95%-CI z p-value
## Fixed effect model -0.0149 [-0.0253; -0.0045] -2.82 0.0048
## Random effects model -0.0109 [-0.0804; 0.0586] -0.31 0.7585
##
## Quantifying heterogeneity:
## tau^2 = 0.0273; H = 5.92 [5.36; 6.54]; I^2 = 97.1% [96.5%; 97.7%]
##
## Test of heterogeneity:
## Q d.f. p-value
## 876.35 25 < 0.0001
##
## Results for subgroups (fixed effect model):
## k 95%-CI Q tau^2 I^2
## paradigm = passive 8 0.0177 [-0.0043; 0.0396] 398.01 0.0642 98.2%
## paradigm = active 18 -0.0242 [-0.0360; -0.0125] 467.45 0.0200 96.4%
##
## Test for subgroup differences (fixed effect model):
## Q d.f. p-value
## Between groups 10.89 1 0.0010
## Within groups 865.46 24 < 0.0001
##
## Results for subgroups (random effects model):
## k 95%-CI Q tau^2 I^2
## paradigm = passive 8 0.1046 [-0.0800; 0.2892] 398.01 0.0642 98.2%
## paradigm = active 18 -0.0608 [-0.1339; 0.0124] 467.45 0.0200 96.4%
##
## Test for subgroup differences (random effects model):
## Q d.f. p-value
## Between groups 2.67 1 0.1026
##
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#no difference in random effects model by paradigm but fixed yes. Hetereogeneity is significantly different so a. fixed not representative and b. need a better model
forest(m, xlim="symmetric", plotwidth=unit(1, "cm"))